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Abstract 

A model for a non-Kerr cylindrical nematic fiber is presented. We 
use the multiple scales method to show the possibility of construct- 
ing different kinds of wavepackets of transverse magnetic (TM) modes 
propagating through the fiber. This procedure allows us to generate 
different hierarchies of nonlinear partial differential equations (PDEs) 
which describe the propagation of optical pulses along the fiber. We go 
beyond the usual weakly nonlinear limit of a Kerr medium and derive a 
complex modified Korteweg-de Vries equation (cmKdV) which governs 
the dynamics for the amplitude of the wavepacket. In this derivation 
the dispersion, self- focussing and diffraction in the nematic are taken 
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into account. It is shown that this cmKdV equation has two-parameter 
families of bright and dark complex solitons. We show analytically that 
under certain conditions, the bright solitons are actually double em- 
bedded solitons. We explain why these solitons do not radiate at all, 
even though their wavenumbers are contained in the linear spectrum 
of the system. We study (numerically and analytically) the stability 
of these solitons. Our results show that these embedded solitons are 
stable solutions, which is an interesting property since in most systems 
the embedded solitons are weakly unstable solutions. Finally, we close 
the paper by making comments on the advantages as well as the limita- 
tions of our approach, and on further generalizations of the model and 
method presented. 

PACS numbers: 42.65.Tg, 61.30G, 77.84.N 
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1 Introduction 



Theoretical studies on the existence of solitons in liquid crystals (LCs) started 
in the late sixties and early seventies PQ-[I], and experimental confirmations 
were reported subsequently [E]-[Bj- In the case of static solitons in LCs, the 
molecular configurations may be obtained from the Lagrange equations derived 
from the Helmholtz free energy, whereas for propagating solitons the continu- 
ous change in these configurations makes it necessary to take into account the 
damping of the molecular motion. For liquid-crystal waveguides, the nonlin- 
earity necessary for the existence of solitons is provided by the coupling with 
the optical field. 

Coupling of the dynamics of the velocity and director fields in LCs to ex- 
ternal optical fields renders the relevant dynamical equations highly nonlinear, 
which makes it possible to have solitary waves of the director field with or 
without involving the fluid motion. Furthermore, the strong coupling of the 
director to light makes any director-wave more easily detectable by optical 
methods than it is in isotropic fluids, where only the flow field is observable. 

Some nonlinear partial differential equations (PDEs) appearing in the 
liquid-crystal theory give rise to exact soliton solutions. These are the Korteweg- 
de Vries (KdV), nonlinear Schrodinger (NLS), and the sine-Gordon (sG) 
equations [9 . The KdV equation describes a medium with weak nonlinearity 
and weak dispersion, whereas the NLS equation describes situations where 
weak nonlinearity and strong dispersion prevail, such as the propagation of 
signals in liquid- crystal optical fibers. 

Passing continuous laser beams through nematic LCs reveals the existence 
of static spatial patterns in cylindrical ^H] and planar 11] geometries. The ba- 
sic physical mechanism which support these time-independent patterns is the 
balance between the nonlinear refraction (self-focussing) and spatial diffrac- 
tion in the nematic. However, when the propagation of wavepackets, rather 
than continuous beams, is considered, a different situation occurs. The en- 
velope of the wavepacket obeys an NLS equation, which takes into account 
self- focussing, dispersion and diffraction in the nematic P^j - [To] . This 
equation has soliton solutions whose speed, time and length scales may be 
estimated by using experimentally measured values of the corresponding pa- 
rameters JHl- However, the usual analysis of this situation is based on the 
assumption that the LC behaves as a Kerr medium and that, consequently, 
strong dispersion and weak nonlinearity, at order 0(q 3 ), with respect the field 
amplitude q, should be taken into account As it will be discussed below, q 
measures the ratio of the electric-field energy density and the elastic-energy 
density of the nematic and it is, therefore, a measure of the coupling between 
the optical field and the fluid. However, although truncating the analysis at 
the 0(q 3 ) order may be a very reasonable assumption for solid-state optical 
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media, the soft nature of the LCs suggest that the neglect of higher-order 
contributions may not necessarily be a good assumption in this case. 

Recently, the formation of spatial solitary waves in nematic LCs with at 
the light-power level of a few milliwatts has attracted a good deal of interest 
|17j . [TS] . [T§] . [2*U] . It has been experimentally shown that the nonlinearity of 
these media can support solitons in LC line waveguides [21], [22]- 

The main purpose of the present work is to develop an approach that al- 
lows to generate PDEs which describe the propagation of optical pulses in 
nematic LC waveguides beyond the weakly nonlinear limit corresponding to 
the Kerr medium. More specifically, we show that to 0(q A ), and assuming 
that attenuation effects are small, the evolution of the amplitude of propagat- 
ing transverse-magnetic transverse-magnetic (TM) modes is governed by an 
equation with a derivative nonlinearity, which is the complex modified KdV 
(cmKdV) equation, 

u z - euut - 7 \u\ 2 u t = 0. (1) 

see Eq. (JH2J) below. 

The paper is organized as follows. In Sec. 2 we introduce a model of a 
cylindrical nematic cell and set up basic coupled equations for the orienta- 
tional and optical fields. We formulate an iterative procedure to expand these 
equations in terms of the coupling parameter q, which leads to a specific hi- 
erarchy of PDEs. Then, in Sec. 3 we derive dynamical equations governing 
the evolution of the amplitude of propagating TM modes up to the order 
0(q A ). Rescaling the equations, we show that the standard NLS equation is 
obtained at order 0(q 3 ), and that the equation corresponding to 0(q i ) is in- 
deed the cmKdV equation (JTJ). In Sec. 4, soliton solutions to this equation are 
studied. In particular, it is shown that the equation has ordinary bright- and 
dark-soliton solutions, and a continuous family of embedded solitons (ESs), 
i.e. solitary waves which exist inside the system's continuous spectrum of lin- 
ear waves |2S1 • In Sec. 5 we discuss why the ESs can exist in Eq. (0) without 
emitting any radiation, even though their wave numbers belong to the linear 
spectrum. In Sec. 6 we study the stability of the ESs. We conclude the pa- 
per in Sec. 7, which summarizes the results and compares them to previously 
published ones. We also point out advantages and limitations of our approach, 
and discuss possible ways to generalize it. 

2 The Model and Basic Equations 

We consider a cylindrical waveguide with an isotropic core of radius a, dielectric 
constant e c and a quiescent nematic LC cladding of radius b. The initial 
orientational state is depicted in Fig. 1, where the director field obeys the 
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following axial strong-anchoring boundary conditions, 



n(r 
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An optical beam is launched into the guide and propagates through the LC. 
If the field is strong enough to exceed the orientation-transition threshold, the 
initial configuration is changed by reorienting the director field. We assume 
that the induced reorientation occurs only occurs in the (r, z) plane, so that 
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where e r and e z are the unit vectors of the cylindric coordinates. 

Although the incident beam is neither planar nor Gaussian, the normal 
modes within the cavity are cylindrical plane waves propagating along the 
z axis. In previous works it has been shown that only the TM modes, with 
nonzero components E r (r, z, t), E z (r, z, t) and H^r, z, t) of the electromagnetic 
field, couple to the reorientation dynamics of the director field [H], ^2], |24j . 
As it can be shown that E r (r, z, t) and E z (r, z, t) may be expressed in terms of 
H^rjZjt), below we only describe the dynamics of the component H^(r, z,t). 
The relevant dynamical equations, which take into account retardation effects, 
are given by Eqs. (8) and (9) of Ref . [21] , namely, 
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with 



S (r*,t) = j j dt j dt'-^j- (V -t'^hh-VxH • (6) 

In these equations, we have used dimensionless variables, £ = z/a, x = r/a, 
= H^/(ceoE ), Ef = S^/E , % = r,z, where E is the amplitude of the 
incident field. The speed of light in vacuum is c = , where /x and 

eo are, respectively, the magnetic permeability and electric permittivity of free 
space. The dielectric anisotropy of the nematic, e a = eii — ej_, is defined in terms 
of the dielectric constant for directions parallel (ey) and perpendicular (e±) to 
the director. As mentioned in Sec. 1, q 2 = eoE^a 2 / K is the dimensionless ratio 
between the electric-field energy density and the elastic-energy density of the 
nematic, where K is its elastic constant in the equal constants approximation. 
Thus, q 2 is a measure of the coupling between the optical field and the LC. We 
stress that, in writing Eqs. (jlj) and (JHJ), large difference between the time scales 
of slow reorientation dynamics and rapid variations of the electromagnetic field 
was explicitly taken into account, and as a consequence the time derivatives 
of 6 were ignored. 

When the coupling between the TM mode H<p(r, z, t) and the reorientation 
field 9(r,z,t) is negligible (q — 0), the propagating modes are represented by 
quasi-planar waves. However, if the nonlinearities in Eq. (j3J) are taken into 
account by considering finite q, they cause space and time variations of the 
field H^(r,z,t), due to generation of higher-order harmonics which feedback 
to the original modes. 

We assume that the interaction between the optical field and the reori- 
entation in the nematic is stronger than in the weakly nonlinear limit (Kerr 
medium) which corresponds to q = 1 12J. Furthermore, in all the analysis 
we neglect all the backflow effects associated with the reorientation or caused 
by external flows [T5] . Thus, we solve the coupled equations ()H) and (J3J by 
assuming the following coupled expansions of 9 and in powers of q, 

9 = 9 {o) + q 2 \A(Z,T)U(x,Lu)\ 2 9 (1 \x) + q 4 \A(E,T)U(x,cu)\ 4 9 (2) (x) + (7) 



H+faCt) = qU^x^o + iX-^Ai-^ + q'uW+q'UW 

+g 4 f/ (4) + g 5 [/ (5) + c.c. + ..., (8) 

where c.c. stands for the complex conjugate. 

The rationale behind this assumption is the following. As indicated in 
Eq. (Jlj), the lowest-order coupling between 9 and occurs at order q 2 , and 
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it is therefore reasonable to expect that higher-order terms will also be even 
m q. The fields 6 {n) with n — 0, 1, 2.. are contributions to 6 at order n , 
which satisfy the same hard-anchoring homeotropic boundary conditions as 
were given above by Eq. (J2J), 9(x = 1) = 9(x = b/a) = 0. As usual, the 
amplitude A(E, T) in Eqs. (JTJ) and (JHJ), which represents an envelope of a 
narrow wavepacket of width A = (u — uj )/uj , whose central frequency is u>q, is 
assumed to be a slowly varying function of the variables H = X( and T = Xt . 
Here A is a small parameter which measures the dispersion of the wavepacket. 
In Eqs. (JTJ and (jHJ), (x.ujo) is the well-known linear solution for which 
is given explicitly by [2H] 



U,A.v.^ () ) J 2 ( \/e c (— ) 2 - (3 2 a 2 ) J ^—exp(-i(3aC - -/ax), (9) 



with 7 = yej| (f3 2 /e± — (u /c) 2 ). Here J\{x) is the Bessel function of order 1, 
and j3 is the propagation constant, which only takes allowed values calculated 
in Ref. [21]. The terms proportional to U^ n \ n = 2,3, ... in (jHJ) are contribu- 
tions to the TM modes from to the higher-order optical harmonics that are 
generated by the nonlinearities in Eqs. (0) and (jU). 

Note, however, that the relation between the parameters q and A is not 
unique. For instance, when the wavepacket is very narrow, this relation is A = q 
and up to 0(q 3 ), the expansion leads to the standard nonlinear Schrodinger 
(NLS) equation for A(E,T) (which corresponds to the Kerr medium) [25]. 
|12j . [24J. Therefore the model may be generalized in various ways. Since q 
and A are small parameters, we assume that A = q a with some positive a. 
Then a = 1/2 represents a wider and a = 2 a narrower wavepacket. Note that 
the presence of higher powers of q implies that these higher-order contributions 
are smaller than the dominant term in (JTJ), which describes a small-amplitude 
narrow wavepacket. 

Inserting the expression (JTJ) into Eq.(J5J) and expanding in powers of q, it 
is straightforward to rewrite Eq.(J5J) as 

L{p, to, x)H$ + q 2 F(H^) + q 4 G(H^) = 0, (10) 

where the linear L operator and nonlinear ones, F and G, are defined, respec- 
tively, as 
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Zero-order solutions for the orientation, with 6^ = 0, and first-order ones, 
which gives rise to d^\x), were found in Ref. [21] by inserting expressions 
(J7J) and (jHJ) into Eqs. (HJ), (fTUj) and solving the resulting equations. This way, 
6 { - l \x) turned out to be 



7re^e||X (a 2 — b 2 ) 1 
+ (6 2 -A 2 )+e^V(l-x 2 )}. (14) 

To study the dynamics beyond the Kerr approximation, we need to calcu- 
late the fourth-order terms in Eq. (|7ji. that is, 8^ 2 \x). To this end, we insert 
Eqs. ©, (0) into Eq. (jlj) and expand the result in powers of g up to the fourth 
order. This leads to 



(4X 2 [6 ± 6||(^f) 2 -/?V]-6 ± ) fl(1)/ 
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After substituting 0^(:r) from Eq. JHJ), this equation takes the form 
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x dx dx 2 

In spite of its apparent complexity, this linear differential equation for 9^ 2 \x) 
can be easily solved by imposing the planar strong-anchoring boundary condi- 
tions for 9, as explained above. The solution can then be written in terms of 
the exponential-integral function, and if the resulting expressions are approx- 
imated by asymptotic expressions for this function, we obtain 



9 {2 \ 



x 



/3eae -(6+«x)7j2[ v / ec (fiM)2_^ 2] 
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While this compact form for 9 iy2 \x) is sufficient for our discussion below, ex- 
pressions for the coefficients A , Ax, A 2 , A 3 , A*, -Bo, B x , C , C\, C 2 , C 3 and 
C4 that appear in Eq. (fT7|) are given in Appendix A. 

To conclude this section, it is relevant to stress that the above derivation 
ignored dissipative loss in the LC medium. In fact, the physical condition 
for the applicability of this assumption is that the propagation distance to 
be passed by excitations (solitons) is essentially smaller than a characteristic 
dissipative-loss length. This condition can be readily met in situations of 
physical relevance. 



3 The Envelope Dynamics 

We now aim to derive an equation for the envelope A(H, T) by dint of the same 
procedure that was used in Ref. ^2] for the weakly nonlinear case. To this 
end, we substitute Eq. (JHJ) into Eq. ([TUj) . and identify the Fourier variables 

i(3a = i^a + q a -^r + + q" a ^~ + , (18) 
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-iu = -iujQ + 1 a gj;, (19) 

where the variables H n , n = 1,2, 3, 4, are related to the spatial scales associated 
with upper harmonics contributions, that is, Z = g raa S n . This substitution 
leads to an equation, 



= L[z{3 a + q a — + q 2a — + q 3a — + q* a — ,-tLU (20) 
d^i d^ 2 o^ 3 d^ A 

+q a -^}H^x, C, t) + q 2 F[H^ (x, (, t)\ + q 4 G[(H^ (x, (, t)\. (21) 

We now fix a = 1, which means selection of the type of the wave packet 
to be considered; choosing a = 2 or a = 1/2 would imply, respectively, a 
narrower or wider packet of the TM modes than for a — 1. In this case, we 
collect contributions to the same power of q, arriving at the expressions 

q 1 : L (i(3 a, -ioJ , %) 17$ (x, oj ) A = 0, (22) 



(F)TJ f) 3 f) \ 



(24) 



q 3 : LU^ = S 2 U^j A - F(U, (x, u> ) A), 

q ■ LU "H^'^'^' 7 

(25) 

\ ouj ox oujox dx J 

where L n , n = 1,2, denotes the derivative of L (i/3 Q a, —iujq) with respect to 
its first or second argument. Clearly, the same procedure can be carried out 
for a = 1/2 or 2. 

Note that Eq. ()22|) is actually the usual dispersion relation L (i(3 a, —iuj ) 
Tibial, Uq) = 0, which confirms our approximation, since H^(x,cuo) already 
satisfies this equation to the first order in q. To simplify Eqs. (|2*3*j) - (|23|) we 
take the first four derivatives of Eq. (|22j) with respect to u. This leads to a set 
of linear inhomogeneous equations for , the existence of solutions to which 
is secured by the so-called alternative Fredholm condition j22j . This condition 
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is fulfilled if LU^ (x, u ) = and if (x, u ) — > as x — > oo. In our case, this 
reads explicitly 

LU {n \ u <t>) = ^ U^LU {n) dx = 0, n = 1, 2, 3, 4. (26) 

By applying the relations (|2l)j) to Eqs. (|2l?jl - (|23|) . substituting the four first 
derivatives of Eq. ()22)1 into them, and collecting terms in front of the same 
power of q, we obtain the following equations for A(S, T) on each of the spatial 
scales S, Si, S 2 , H 3 , S 4 , for the successive orders in q, 



9 : 5^-5S?ffl5- 9T = °- (29) 

Here, dimensionless coefficients n 2 and n 3 are defined as follows: 
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The coefficient Hi is related with the nonlinear diffraction index n 2 through the 
expression n 2 = Kn^j^a 2 ■ Similarly, we define a nonlinear diffraction index 
at the next order beyond the Kerr approximation by 713 = uo^Kn^j e^a 2 ] it is 
proportional to the coefficient in front of the nonlinear term in Eq. (|29j) . Note 
that Eq. ()27j) simply describes a wavepacket in the linear medium, while Eq. 
(jZSj) is the well-known NLS equation which gives rise to robust soliton pulses. 
The equations corresponding to the orders q 2 and q 3 are well-known ones, and 
they have also been derived and analyzed in Ref. [24J. In the next section we 
focus on Eq. (|2^|) . which was derived at order q 4 . 
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4 Double Embedded Solitons 



Equation (J29J) may be rewritten in a rescaled form by introducing the dimen- 
sionless variables u = A/A , £ = H 4 /Z 4, and r = T/T 4, where Z Qi and T 4 
are space and time scales, and A is the initial amplitude of the optical pulse. 
In terms of these variables, Eq. (|2H|) becomes 

du d 3 u 2 du 

d~r £ d^~ llul ^ =0, (32) 

where we have defined the dimensionless coefficients e and 7 as 

1 Z 04 d 3 (3 



6 T 04 3 rfw 3 ' 



(33) 



7 = /3n 3 A 2 ^i. (34) 

J 04 

In what follows below, we will consider Eq. ()32|) in the form of Eq. i.e., 
with £ and r replaced by 2 and t. 

Equation (|3*2|). or equivalently Eq. (JTJ), reduces to the real modified Ko- 
rteweg de Vries (mKdV) equation when we restrict u(z,t) to be real, hence 
all the real solutions of the mKdV equation, including iV-soliton ones, are 
also solutions of On the other hand, Eq. (jlj also has complex solutions 
which include, as it will be discussed below, two-parameter families of bright 
and dark complex solitons. Actually, the existence of these complex solutions 
of Eq. (JTJ) was pointed out by Ablowitz and Segur as early as 1981 [2H]- The 
precise form of the bright solitons in the particular case when e = 67 was 
presented recently by Karpman et al. [3Q|. 

In the general case the bright-soliton solutions to Eq. (0) may be found by 
substituting a straightforward trial function in this equation, 

\ J(qz+rt) 



u(z, t)=A sech [—^- j e t(qz+rt) . (35) 
This substitution shows that (|33j) is indeed a solution of (JTJ), provided that 

A 2 w 2 = -, (36) 

7 

a = 3er 2 - ^ 7 A 2 , (37) 


q = \pA 2 r - er 3 . (38) 
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Condition (|3*S)l implies that the bright soliton solution ([3*3|l only exists for 
67 > 0, which implies that, in the opposite case, the nonlinearity and linear 
dispersion cannot be in balance. Moreover, since we have five free parameters 
in Eq. and only three conditions - (|5Sjl. these expressions define 

a two-parameter family of bright soliton solutions of Eq. (^J), so that the 
following pairs of the parameters can be chosen arbitrarily: (A,r), (w,r), 
(A,q), or (w,q). The family includes, as particular cases, the real one-soliton 
solutions of the mKdV equation, which are obtained when r = 0. 

In a similar way, dark solitons of Eq. ((H) can be found by substituting the 
trial function 

u(z,t) = A d tanh ( ^* +r " t \ (39) 
\ w d J 

This substitution shows that this ansatz solves Eq. if the following condi- 
tions are satisfied 

4^ = --, (40) 

7 

a d = 3er d 2 -^Al (41) 

q d = 1 Ajr d -erl, (42) 

which are similar to the conditions (|36|) - (|38|) for the bright solitons. As in 
the bright case, the conditions (|3U|) - permit us to choose freely any of 
the following pairs of parameters: (A, r), (w,r), (A,q), or (w,q). Thus, Eqs. 
(|39|) - (j42|) define a two-parameter family of dark-soliton solutions of Eq. (fT|). 
Eq. (|4()jl showing that this family only exists if £7 < 0, i.e., exactly in the case 
opposite to that in which bright solitons are found. 

Out of the two families of the above soliton solutions (bright and dark) 
of Eq. (Q), the bright family is the most interesting one. In spite of their 
similarity to ordinary bright solitons, the bright soliton solutions of Eq. (|3*2*j) 
feature a special property which distinguishes them from ordinary solitary 
waves, namely, they are double- embedded solitons. The concept of embedded 
solitons (ESs) was formulated, in a general form, in Ref. It refers to 

solitary waves which do not emit radiation, in spite of the fact that the soliton's 
wavenumber (spatial frequency) is embedded in the system's linear spectrum. 
Still earlier, solitons of this type were found in particular models |17j . for 
instance, in a generalized NLS equation involving a quintic nonlinear term 
[3*T] . Recently, more systems supporting ESs have been found [32] - |HI]- To 
the best of our knowledge, the existence of ESs has not been reported before 
in models of LC media. 
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So far, the embedded solitons were classified in two groups, namely, those 
which obey NLS-like equations (or systems thereof), and those which are 
governed by KdV-\ike equations. In the former case, an ES has its wavenumber 
embedded in the range of wavenumbers permitted to linear waves (as it was 
already mentioned above). In the latter case, the velocity of an ES is found 
in the range of phase velocities of linear waves. There are, accordingly, two 
different ways to decide whether a solitary-wave solution to a nonlinear PDE 
system is embedded, viz., the wavenumber (WN) and velocity (VE) criteria. 

In Ref. jHDl it was pointed out that Eq. ((H) is a particular case of a more 
general iVLS'-like equation possessing ESs. For this reason, and also in view 
of the significance of Eq. for physical applications, it is interesting to 
determine if the bright-soliton solutions of Eq. may be ESs. It should 
be noted that Eq. (JJJ) may be regarded as both a KdV-\ike equation, due 
to its similarity to the mKdV one, and an NLS-Wke equation, because, in 
the context of wave propagation in LCs, Eq. (0) in its complex form plays 
a role similar to that of the NLS equation, i.e., the one governing evolution 
of a slowly varying envelope of a rapidly oscillating wave. Therefore, it may 
be possible to apply both criteria, WN and VE ones, to decide if the soliton 
solutions of Eq. (|JJ) are ESs. 

First, we apply the WN criterion. To this end, we must determine if the 
wavenumber of the solution ()35|) is contained within the range of the wavenum- 
bers allowed to linear waves. To identify the intrinsic wavenumber of the so- 
lution, we must transform it into the reference frame moving along the time 
axis with the reciprocal velocity a, see Eq. (135)1 . The transformation adds a 
Doppler term to the soliton's internal spatial frequency (wavenumber), making 
it equal to q + ar. On the other hand, plane-wave solutions to the linearized 
version of Eq. (JJJ) in the same reference frame can be sought for as 

u(z, t) = exp i [kz — uj{t — az)} , (43) 

which leads to the following dispersion relation 

k{u) = eu 3 - au. (44) 

Since the range in which the function (|44jl takes its values covers all the real 
numbers, including the soliton's wavenumber q + ar, all the soliton solutions 
to Eq. (JJJ), given by Eqs. (|33 jl - (|S5jl . are classified as ESs as per the WN 
criterion. 

Now, we address the question whether these solitons are also embedded 
according to the VE criterion. As the evolution variable in Eq. (JJJ is the 
distance z, rather than the time t, it is the reciprocal velocity which determines 
if the moving solutions are embedded according to the VE criterion. Thus, 
we should find out if the reciprocal velocity of the soliton (|35jl . given by the 
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parameter a, is contained within the range of the reciprocal velocities permitted 
to linear waves. The dispersion relation (|44|) implies that the reciprocal phase 
velocities of the linear waves (in the reference frame moving along with the 
soliton) are given by 

- = -a + eu 2 , (45) 
us 

while the reciprocal velocity of the soliton proper is, obviously, zero in the 
same reference frame. Obviously, the expression ()45|) takes the value zero if 
ae is positive, hence the soliton solutions given by Eqs. (J35|) - (J35j) are ESs 
according to the VE criterion provided that ae > 0. As these solitons are 
also embedded according to the WN criterion, we call them double- embedded 
solitons. On the other hand, when ae < 0, the soliton solutions of Eq. are 
only embedded with respect to the WN criterion, but not as per the VE one, 
therefore in this case we apply the term single- embedded solitons. 

5 Radiation Inhibition and Continuity of the 
Embedded Solitons 

As in any other system with ESs, the fact that the solitons do not emit 
radiation despite being embedded in the linear spectrum should be explained. 
Since the wavenumber q + ar of the soliton solution (jH5)l is contained in the 
linear spectrum defined by the dispersion relation (j44|) . a resonance of the 
soliton is expected with the linear waves whose frequencies satisfy the condition 

q + ar = ecu 3 — acu. (46) 

Moreover, when ae > the soliton's reciprocal velocity a coincides with the 
reciprocal phase velocities (eu 2 ) of two linear waves whose frequencies satisfy 
the condition 

a = ecu 2 , (47) 

consequently one could also expect the soliton to resonate with these waves. 
Different explanations for the absence of resonant radiation in other systems 
which support ESs were proposed [33], jH]. However, an explanation for the 
radiationless character of the ESs in Eq. has not been presented. 

Another unexpected property of the same ESs in Eq. is the fact that 
they exist in a continuous family. In most cases, ESs are isolated solutions; 
usually they do not appear in families, although examples of continuous fam- 
ilies of ESs are known too, for instance, in a fifth-order KdV equation 35 . 
It is also necessary to explain why Eq. (JIJ) has a two-parameter family of the 
ES solutions. 
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As we show below, the radiationless character of the ESs in Eq. (0) is the 
consequence of a special balance between the linear and the nonlinear terms 
of this equation. To understand how these terms interact, it will be helpful to 
separate their effects by considering the following linear driven equation, 

du <9 3 w 2 du 

d- z - £ w~ llu ° l ^r = ' (48) 

where the source is built of a solution Uo(z,t) to Eq. (pQ). It is clear that the 
same function uq is also a solution to Eq. (flfijl . 

We now define the double Fourier transform of u(z,t), 

oo oo 

u(k,Lu) = ^-J J ' u(z,t)e- i{kz - ut) dzdt, (49) 

— oo — oo 

and Fourier transform Eq. (J48j) . to obtain 

u(k,u)=i ^^3 , (50) 

where 

F (z,t)= 1 \u \ 2 ^ (51) 

is the source in Eq. pSjl. 

To understand the mechanism of the cancellation of the emission of ra- 
diation, we can temporarily take, instead of the exact source (|51J1 . a model 
example with 

F (z,t) = Asech y-^-J e i{qz+rt) . (52) 

In this case, the calculation of the Fourier transform of F Q and substitution in 
Eq. (J5UJ) yield a result 

, irwAsech. \%w(r + tv)l { A 2 ^r w 2 A 2 'jr 3 A 2 ^uj 

U{k,(jJ) = < h 



(r + u)a — q + ecu 3 [ 3 3 6 

w 2 A 2 ^r 2 u w 2 A 2 ^u' 



+ ^—>S{[(r + u)a + q]-k}. (53) 

At first sight, this expression seems to imply that a resonance with the ra- 
diation waves should occur for frequencies at which the denominator, which 
is a third-order polynomial in u, vanishes, which is actually tantamount to 
Eq.(JHI]). Moreover, if q — r = 0, the same argument shows that a resonance 
at the frequencies defined by Eq. (|37j) should be expected. Observe, however, 
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that the numerator on the right-hand side of (|53|) also contains a third-order 
polynomial in uj. Consequently, if the two polynomials happen to coincide, 
they will cancel each other, which also implies the cancellation of the resonant 
generation of the radiation modes. Equating the coefficients in front of powers 
of uj in the two polynomials in ([53)1 . we obtain three equations which, after 
some manipulations, take the precise forms of Eqs. ([3TI|) - (pIHj) . Thus, these 
three equations are the necessary and sufficient conditions for the mutual can- 
cellation of the two polynomials in Eq. ([53)1 . This explains why the forcing 
term F (z,t) of the form (}52"j) does not generate any radiation, provided that 
the parameters A, a, w, q and r satisfy Eqs. (jHEJ) - ([HE]) - Furthermore, observe 
that the polynomial that appears in the numerator of the expression (153)1 
contains the nonlinear coefficient 7, while the polynomial in the denominator 
contains the dispersion coefficient e. Consequently, the cancellation between 
these two polynomials is a result of the balance between the nonlinearity and 
dispersion in Eq. Finally, if the forcing term in Eq. ([57I|) is taken in the 
exact form ([51)1 . rather than in the simplified form of Eq. ([52j) . it is easy to 
check that the prediction for the cancellation of the radiation emission will be 
the same. 

In the case of the full equation (JT)), the same cancellation argument explains 
why an initial condition of the form 



does not radiate at the frequencies defined by Eqs. (}46|) and ([4?]) if A and w 
do not satisfy Eq. ()36|) . On the other hand, if A and w this condition, we 
expect the resonances to occur. In the next section, we will verify numerically 
that this indeed the case. 

To close this section, it is relevant to stress that the cancellation of the 
two polynomials in Eq. ([H3*|) imposes only three conditions, while the so- 
lution ([33|l involves five parameters. Therefore, the cancellation conditions 
do not uniquely determine the soliton parameters, which explains why the 
soliton solution ([33|) - ([3*K[) involves two arbitrary parameters, thus defining a 
two-parameter continuous family of the ESs. 

6 Stability of the Embedded Solitons 

In this section we will study the stability of the bright-soliton solutions of Eq. 
(£Q). As it was explained in Sec. 4, these solitons may be either single- embedded 
or double- embedded, depending on the sign of the parameter combination ae. 




(54) 
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In the following we will separately consider the cases of positive and negative 
as. 

We begin by considering a single-embedded soliton of Eq. (0) , setting e — 1 
and 7 = 6 [these values were chosen as they correspond to those at which the 
related Hirota equation [12], which is connected to Eq. by the Galilean 
transform j^Jjj, is an exactly integrable one [IB]]. We start with the following 
values of the soliton parameters 



0.790 (55) 




(56) 
(57) 



-1/2 (5* 



11 , s 

q s = — =. (59) 



6V24 

These values satisfy the conditions (}3l)j) - (J3~%j) . and therefore they characterize 
an exact bright soliton of the form ()35|). Since a s e < 0, this soliton is a single- 
embedded one (i.e., it is embedded solely according to the WN criterion). 
To test stability of this soliton, we consider an initial condition of the form 

u(z = 0, t) = A sech ( — ) exp (ir t) , (60) 



where w = w s and r = r s , but A is slightly different from A s . If we give A 
a value 0.815, which is larger than A s , the numerical solution of Eq. shows 
that the pulse moves to the right along the temporal axis with a reciprocal 
velocity equal to —0.58, which is slightly lower than a s , and the pulse's ampli- 
tude evolves as shown in the upper curve of Fig. 2. The observation that the 
reciprocal velocity of the perturbed pulse is lower than a s is consistent with 
Eq. (jSTj), which indicates that a should decrease if A is increased. Figure 2 
(the upper curve) shows that the pulse's amplitude stabilizes and approaches 
an equilibrium value close to A = 0.84. The temporal profile of the pulse 
at z = 50 is displayed in Fig. 3. It shows a small-amplitude radiation wave 
emitted by the trailing edge of the pulse. The frequency composition of this 
tiny radiation wave can be determined by calculating the Fourier transform 
(FT) of the radiation contained in the interval 40 < t < 168. The insert in 
Fig. 3 shows the power spectrum (i.e., the square of the FT amplitude) of 
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this radiation, which contains two peaks located at the frequencies v\ = —0.10 
and z/ 2 = 0.12. These peaks are close to the resonant frequencies {y = ±0.07) 
predicted by the resonance condition 



These radiation peaks imply that the perturbed pulse emits radiation according 
to the way the soliton's wavenumber is embedded in the spectrum of linear 
waves. 

If we now consider an initial condition of the form (|60|). with w a = w s , 
r = r s and A = 0.765 < A s , the behavior of the perturbed pulse is similar. 
In this case the amplitude evolves as shown in the lower curve of Fig. 2 where 
we can see that the pulse's amplitude approaches an equilibrium value close 
to A = 0.74. The reciprocal velocity of the perturbed pulse is —0.42, which 
is slightly higher than a s . This change is consistent with Eq. pTj) . which 
indicates that a should increase if A is diminished. 

The two curves shown in Fig. 2 demonstrate that the single-embedded 
soliton solutions of Eq. (0) are stable. This is an interesting result, since 
usually ESs display a weak (nonlinear) one-sided instability 34 a - In fact, the 
complete stability of the ESs in Eq. (0) may be expected, due to the fact that 
in this case we are dealing with a continuous two-parameter family of the ESs, 
while in most other systems ESs are isolated solutions, which explains their 
nonlinear instability. 

Figure 2 also shows that if the amplitude of one of the single-embedded 
solitons of Eq. ((TJ) is slightly increased, the perturbed soliton stabilizes itself 
at an even higher amplitude. On the contrary, if the soliton's amplitude is 
slightly decreased, the perturbed soliton stabilizes at a still lower amplitude. 
This behavior can be better understood if we analyze the evolution of the 
perturbed solitons of Eq. (0) by means of the averaged variational technique 
introduced by Anderson jH], which is one of the approximately analytical 
methods used successfully in nonlinear optics [IS] |5T] . see also a recent 
review 

In order to apply the variational technique, we start with the ansatz of the 
ordinary form, 



q s + a s r s = euo 3 - a s uj 
[cf. Eq. (j4fij) ] and the partial resonance condition 00] 

- (q s + a s r s ) = eu 3 - a s u. 



(61) 



(62) 




(63) 



L = i (u z u* - u* z u) + ie {uu* ttt - u*u m ) + tt [' 



(u*) 2 uu t ] 



(64) 
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and integrating over time, we calculate the averaged (effective) Lagrangian 

/■oo 

C = Ldt. (65) 



The following Euler-Lagrange equations can be easily derived from C, 

-8AWQ' - - 8sAR 3 W + — jA 3 RW = h{P,P',R'), (66) 

vv 3 



AA 2 Q> + - AeA 2 R 3 + i^R = / a (P, P', A',V',W', R'), (67) 



f 3 (P,P',A',V',W',R') = 0, (68) 

A 2 W = A 2 {0)W{0), (69) 

iOp /1 2 4 

12eA 2 J R 2 ^ + -7A 4 ^ = /4(P,P , ,^ , ,^',W r/ , J R , ) ) (70) 



" " " 3 



24e^ 2 J RViy+ - 7 AV^ = f^P.P'.A'.V, W',R'), (71) 

It 3 

where the primes stand for the ^-derivatives, and the expressions f n (P, P ', A', V ' 
,W',R') are nonlinear functions of their arguments. Their explicit forms are 
not given, as they will not be needed in what follows. 

We now resort to search for fixed points of Eqs. (JBTIj) - (|71|). which are 
stationary solutions of the form 

A' = W' = R' = P' = P = Q, (72) 

Q ' = const = q, (73) 

V' = const = o. (74) 

When we insert these conditions into Eqs. (jo^j) - (|7T|) . we find that f n — (for 
n = 1, 5) and consequently, the following relations are obtained 

A 2 ^ 2 = i^£ (75) 

7 
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a = 3eR 2 - -7A 2 , 
6 



q = l - 1 A 2 R-eR\ 



(76) 
(77) 



A 2 W = const = A 2 (0) W{0). (78) 

Equation (J75J) is the variational counterpart of Eq. (j3Tlj) . and the expressions 
for a and q coincide exactly with those in (|37J) and ()38|). On the other hand, 
Eq. (J7HJ) applies not only to stationary solutions, but to general dynamical 
equations as well, with variable A(z) and W(z), as it expresses the variational 
version of the exact conservation law (which is simply the energy conservation 
in the case of nonlinear optics |S2])- 

Equation (fT8|) is plotted by thin curves, corresponding to two different 
initial conditions, in Fig. 4. This figure also shows plots (the bold curve) of 
Eq. (|7H|) . corresponding to e — 1 and 7 = 6. This diagram helps to understand 
why the soliton [characterized by the parameters - (|59jl and marked by 
point E on the bold curve in Fig. 4], if perturbed by increasing or decreasing 
its initial amplitude, stabilizes itself, as was observed in Fig. 2. 

We take, as the initial perturbed soliton, the one corresponding to point 1 
in Fig. 4. It has the same width as the unperturbed soliton at point E, but a 
larger amplitude, 

A x = 0.815 > 0.790 = A E . (79) 

According to Eq. (J7HJ), the perturbed pulse must evolve sliding along the thin 
curve passing through point 1. The thin curve intersects the equilibrium bold 
curve at point 2, which is therefore a fixed point. Within the framework of the 
variational approximation proper, the trajectory may perform some oscillations 
in a vicinity of this fixed point; however, if effective loss due to the emission of 
small amounts of radiation by the perturbed soliton (which was observed above 
in direct simulations) is taken into regard, the trajectory will be attracted to 
the fixed point, and will eventually end up being trapped at this point, thus 
implying the stabilization of the soliton very close to point 2, which has the 
value of the amplitude 0.840. 

Similarly, starting at the initial condition corresponding to point 3, the 
soliton will slide along the thin line until it gets stuck at the stable fixed 
point 4. As the amplitude corresponding to point 4 is 0.740, the origin of 
the stabilization process observed in direct simulations displayed in the lower 
curve of Fig. 2 is now clear. 

So far we considered relaxation of perturbed single-embedded soliton. Now 
we proceed to the stability of double-embedded ones. To this end, we set 
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e = 7 = 



1, and choose the soliton parameters 

A d=\ /I ~ 0.790, 



(80) 



w d = 




(81) 



r d = 1/4, 



(82) 



a d = 1/12 » 0.08, 



(83) 



= 1/16. 



(84) 



These values satisfy the conditions (j3Tlj) - (J3SJ), therefore they define an exact 
bright soliton of the form (J33j). Since a d e > 0, this soliton is a double-embedded 
one. We perturb it by taking an initial condition of the form (joTIj) with wq = w d , 
r = r d , and A = 0.815 > A d . 

The numerical solution of Eq. (£[]) corresponding to this initial condition 
shows that the perturbed pulse moves along the temporal axis with a reciprocal 
velocity equal to 0.07, which is slightly lower than a d (this lower velocity is 
consistent with Eq. (JH7J)). Simultaneously, the pulse's amplitude oscillates 
as shown in the upper curve of Fig. 5. This figure again shows a trend 
of the perturbed pulse to stabilize. However, in this case (with the double- 
embedded soliton) the stabilization process is slower, and it is necessary to pass 
a greater distance (along the z axis) to observe the damping of the amplitude 
oscillations. The upper curve of Fig. 5 shows that the pulse's amplitude 
eventually approaches an equilibrium value close to 0.84. 

The trailing edge of the perturbed double-embedded soliton emits a tiny 
radiation wavetrain whose frequency components can be determined by cal- 
culating the FT of the radiation contained in the interval 45 < t < 109 (for 
z = 200). The spectrum obtained in this way is shown in Fig. 6. In this figure 
two peaks are seen. The bigger one corresponds to the frequency v = —0.046 
{oj = —0.289), which corresponds to the negative solution of Eq. ( |4T|). and 
therefore it is a consequence of the resonance of the perturbed soliton with 
a linear wave whose phase velocity is equal to the soliton's velocity. On the 
other hand, the smaller radiation peak is located at u = 0.078 {u = 0.490), 
which is very near to the only real root (u = 1/2) of the resonance condition 



Therefore, the latter peak is due to the fact that the soliton's wavenumber 
Qd + o,d fd is contained in the range of wavenumbers permitted to linear waves. 



Qd + a d r d = eu 3 - a d u. 



(85) 
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As the larger radiation peak (the one at v = —0.046) exists due to the 
fact that the soliton is embedded according to the VE criterion, one could 
assume that in this case (i.e., when a double-embedded soliton is perturbed) 
the radiation emitted by the pulse is mainly due to the VE embedding of the 
soliton. However, such a conclusion would be wrong. The left radiation peak 
in Fig. 6 actually has a larger amplitude because the FT of the complete 
solution is slightly shifted to the left (as a consequence of r being positive), 
and it is this shift which enhances the left radiation peak. 

To verify the latter point, one can consider a slightly different initial con- 
dition, characterized by the parameters A = 0.815, wq = Wd = ^/48/5, and 
ro = — Td = —1/4. As in this case r$ is negative, the FT of the complete solu- 
tion will be shifted to the right, and this shift will enhance the right radiation 
peak. At the insert in Fig. 6 we show the spectrum of the radiation emitted 
in this case by the perturbed double-embedded soliton (for z = 100). As ex- 
pected, in this case the radiation peak due to the wavenumber embedding of 
the soliton (i.e., the right peak) is much higher than the one existing due to 
the velocity embedding (the very small peak on the left). We thus conclude 
that both embeddings, WN and VE, are important to explain the emission of 
radiation by perturbed double-embedded solitons. 

If we now consider an initial pulse of the form (JoTij) with Wq = Wd, r = r& 
and A = 0.765 < A d , the numerical solution of Eq. (JTJ) shows that the pulse's 
amplitude again performs a damped oscillatory behavior, as shown in the lower 
curve of Fig. 5. In this case, the pulse's amplitude approaches an equilibrium 
value close to 0.74. 

7 Concluding Remarks 

In this work, using the multiple scales method, we have derived a model for the 
propagation of a wavepacket of TM modes along a cylindrical liquid-crystal 
waveguide beyond the usual weakly nonlinear limit of the Kerr medium. In this 
case, the amplitude of the wavepacket obeys a nonlinear equation, ((H) or (|3^j) . 
which exhibits a derivative nonlinearity. This complex modified KdV equation 
gives rise to the two-parameter families of bright, Eqs. (j3*B]) - (jUBJ), and dark, 
Eqs. (139|) - (j42j), solitons. The bright-soliton solutions of Eq. (JTJ are embedded 
solitons (ESs) (or sometimes double-embedded ones), i.e., they do not emit 
any radiation, in spite of the fact that their wavenumbers (and sometimes their 
velocities too) fall into the linear spectrum of the system. We have shown that 
the physical nature of the existence of the ESs inside the continuous spectrum 
is the balance between the dispersion and nonlinearity in Eq. . Moreover, it 
was concluded that these ESs are completely stable solutions, while, in most 
previously considered models, they are weakly unstable. It was observed that 
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perturbed single-embedded solitons relax to a new equilibrium state faster than 
double-embedded ones. 

The coupled expansions for 9 and in powers of q, that were introduced 
in Sec. 2, can be extended to higher orders. This leads to nonlinear equations 
with the quintic i.e., 0(q 5 ), nonlinearity. Investigation of the corresponding 
model is currently in progress. Also, as discussed in Sec. 3, up to the order 
0(q A ) considered here, the same procedure to construct narrower (a = 2) or 
wider (a = 1/2) wavepackets of TM modes can also be carried out. 

Another possible generalization of our model, not dealt with here, is a 
possibility to take into account hydrodynamic flows beyond the Kerr-medium 
approximation, that will inevitably couple to the reorientation dynamics of 
the liquid crystal. Actually, the inclusion of the flow is unavoidable owing 
to the fluid nature of the system. However, the consideration of the hydro- 
dynamical part of the system substantially complicates the problem. Some 
effects produced by this generalization were considered, at the level of the 
NLS approximation, i.e., at order 0(q 3 ), in Ref. (THj. 
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Figure Captions 

Fig. 1. Schematic of a laser beam propagating through the nematic liquid- 
crystal cylindrical guide. Transverse-magnetic (TM) modes are shown explic- 
itly. 

Fig. 2 Evolution of the amplitude of two perturbed single-embedded soli- 
tons of Eq. ((TJ) (with e = 1 and 7 = 6). The upper curve corresponds to 
the initial condition (j6T?j) with A = 0.815 > A s , w = w s and r = r s , where 
A s , w s and r s are the values (|B3j) - (J57j) . The lower curve corresponds to a 
similar initial condition with A = 0.765 < A s [A(z) and z are dimensionless 
quantities]. 

Fig. 3. Temporal profile (at z = 50) of the perturbed single-embedded 
soliton of Eq. ((TJ) whose amplitude (as a function of z is shown in the upper 
curve of Fig. 2. The spectrum of this profile is shown in the insert [u and t 
are dimensionless quantities]. 

Fig. 4. The bold curve passing through the point E = (Ae,We) = 
(0.790, 2.192) is the plot of Eq. JZ5J) with e = 1 and 7 = 6. The thin 
line passing through point 1 plots Eq. (fTHj) with A(0) = 0.815 > Ae and 
W(0) = We- The thin line passing through point 3 is also a plot of Eq. (fTHj) . 
with A(0) = 0.765 < A E and W(0) = W E - 
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Fig. 5. Evolution of the amplitude of two perturbed double-embedded 
solitons of Eq. ((H) (with £ = 7 = 1). The upper curve corresponds to the 
initial condition (joTJj) with A = 0.815 > A d , w = w d , and r = r d [where A d , 
w d and r d are given by Eqs. (jHUj) - ([52])]. and the lower curve corresponds to a 
similar initial condition with A = 0.765 < A d [A(z) and z are dimensionless 
quantities]. 

Fig. 6. Spectrum (obtained at z = 200) of the radiation emitted by the 
perturbed double-embedded soliton whose amplitude (as a function of z) is 
shown in the upper curve of Fig. 5. The insert shows the spectrum (obtained 
at z = 100) of the radiation emited when the sign of r d is reversed (i.e., when 
r d = -1/4). 
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Appendix A 

Expressions for the coefficients A , A%, A 2 , A 3 , A4, B , Bi, C , G\ y C 2 , C 3 
and C4, which appear in Eq. (JTTj): 

A {a,b;j,e ± ) = Aa 3 b{a + b)^e ± , (86) 

v4i(a,6;/5,7,/i,e a ,e_L,/c ) = a 2 {-24ab(a + b)[(3 2 e a - e ± (e a + e ± )/j,k%) + 

7 [16a(3a - b)b 2 (3 2 e a - 2(a 2 + ao + 86 2 )e ± + 
16a6 2 (— 3a + o)e_i_(e a + e±)fik 2 ] — 
ab{a + 6)e ± 7 2 }, (87) 



A 2 (a, 6; /5, 7, fjt, e a , e ± , k ) = -12ab(a + b)j[-e± + 4a 2 (/3 2 e a - e ± (e a + e±)fikl)], 

(88) 

A 3 (a,6;/5,7,/i,e a ,e_L,fc ) = a{246(a + 6)(/5 2 e a - e ± (e a + e^)^) + 

2 7 [32a 2 6/3 2 e a + 8ab 2 p 2 e a + 8b 3 p 2 e a + 
ae ± - 7be ± - 86(4a 2 + ab + b 2 )e ± (e a + e±)[ik%] + 
7 2 6(a + 6)e ± }, (89) 



A 4 (a, 6; /?, 7, /i, e a , e ± , k ) = 16a6(a + 6) 7 [-/3 2 e a + e±(e a + e±)fik%], (90) 
B (a, b; 7, e ± ) = -2a6(a - 6) (a + 6) 2 7 e_L, (91) 



B x (a, 6; /?, 7, //, e a , e ± , fc„) = ab(a -b){a + b) 2 [2Af3 2 e a - 24e ± (e a + ejj/ifcjj + 7 2 e±] , 

(92) 

C (a, 6; 7, ei) = 4a6 3 (a + 6)7e ± , (93) 

Ci(a,6;/5,7,/i,e a ,e_L, fc ) = o 2 {24ao(a + o)[/5 2 e a - e ± (e a + e±)[ikl] + 

7[-16a 2 6(a - 3b)(3 2 e a - 2(8a 2 + 3ab + 3o 2 )e ± 

+16a 2 o(a — 3o)e_i_(e a + 6j_)/x/cq 

+ao(a + 6)7e ± }, (94) 
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C 2 (a,b;f3,-f,n,e a , ej_, k ) = -V2ab(a + b)j[-e± + Ab 2 ((3 2 e a - e±(e a + e±)iik 2 l )], 

(95) 

C 3 (a,6;/3,7,/i, e a ,e±, k ) = b{2Aa(a + 6)[-/3 2 e a + ej_(e a + e ± )/ik 2 l ] + 

2-f{8a 3 (3 2 e a + 8a 2 b(3 2 e a + 32ab 2 (3 2 e a - 5ae ± 
+3be ± - 8a(a 2 + ab + Ab 2 )e ± (e a + e±)(j,kl) 
-a(a + o)e ±7 2 }, (96) 



C 4 (a, b; (3, 7, /i, e a , e ± , k ) = 16ab(a + b)-f[-f3 2 e a + e ± (e a + e ± )/iA; 2 ] . (97) 
An explicit form of 6^ 2 \x): 



e^ 2 \x) = 



24x 2 a 2 b(a + b) 2 (b 



a)7r 3 7 2 ej]_e 2 



[_ a 2 e (i+x)a 7 (1 _ x) ^2 {24baxf3 2 e a 



ujQa 2 



-48bx[3 2 >ye a (b - ax) + 2(xa - 2b)<ye ± + baxe±>y 2 + 246a;a(— ) 

r 

(-1 + 267 - 2xa7e_Le||)} + bxa {8bax(3 2 e a (3 + 267 - 2xa7) 



,u a. 



+2e±j(xa — 2b + r ybax/2) — 8bxa(3 + 2bq — 2xa r y)ej_e\\ 
+ab(a -b)(a + o)V (6+ah [-2 7 e ± + xa (24(3 2 e a 
+e ± ( 7 2 - 24 e ,|(^) 2 ))] + b(xa - b)e {b+ax) ~< (Qba X1 e ± (b + xa) 

+lQa 4 x(b + xa) 7 [/5 2 e a - e,, (^) 2 ] ) + a 2 (-8e a xa(3 2 [3b + 6o 2 7 

+xa(3 + 2xa^)] — 7ej_ [4(6 — 3xa) + xa^{b + xa)] + 

+8rm [36 + 66 2 7 + xa{3 + 2:ra7)] ej_e|| (— — ) 2 j + a (— 10x 2 a 2 7ej_ 

+46 2 [-6za/3 2 e + 12a; 2 a 2 /3 2 7e a - 7e ± - 7 2 xae ± /4 

— 6xa(2xa / -f — l)ej_e|| ( — ^ — ) 2 + Abxa [— 6xa(3 2 e a — 4x 2 a 2 (3 2/ -fe a 



+7e_i_/2 — 7 2 xaej_/4 + 2xa(2xa7 + 3)e_i_e||( 



uj a 2 



(98) 
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